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ABSTRACT 


Gobiids are hyperdiverse compared with other teleost groups, with about 2,000 species occurring in marine, 
freshwater, and blackish habitats, and they show a remarkable variety of morphologies and ecology. Testing 
the effectiveness of DNA barcodes on species that have emerged as a result of radiation remains a major chal- 
lenge in evolutionary biology. Here, we used the cytochrome c oxidase subunit 1 (COI) sequences from 144 
species of gobies and related species to evaluate the performance of distance-based DNA barcoding and to 
conduct a phylogenetic analysis. The average intra-genus genetic distance was considerably higher than that 
obtained in previous studies. Additionally, the interspecific divergence at higher taxonomic levels was not 
significantly different from that at the intragenus level, suggesting that congeneric gobies possess substantial 
interspecific sequence divergence in their COI gene. However, levels of intragenus divergence varied greatly 
among genera, and we do not provide sufficient evidence for using COI for cryptic species delimitation. Signi- 
ficantly more nucleotide changes were observed at the third codon position than that at the first and the second 
codons, revealing that extensive variation in COI reflects synonymous changes and little protein level variation. 
Despite clear signatures in several genera, the COI sequences did resolve genealogical relationships in the 
phylogenetic analysis well. Our results support the validity of COI barcoding for gobiid species identification, 
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but the utilization of more gene regions will assist to offer a more robust gobiid species phylogeny. 
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INTRODUCTION 


Gobies (family Gobiidae; suborder Gobioidei) are incompa- 
rable among vertebrates in their capacity to adapt and diver- 
sify, which has led to adaptive radiation and rapid speciation 
(Zander, 2011). Gobiid fishes are hyperdiverse compared 
with other teleost groups, with approximately 2,000 species 
in 210 genera occurring in marine, freshwater, and blackish 
habitats. These fish show remarkable morphological and eco- 
logical variety (Nelson, 2006; Zander, 2011). Gobiid fishes 
are globally distributed (Nelson, 2006) and frequently repre- 
sent a dominant component of coral reefs and coastal fish 
communities throughout much of their range, accounting for 
> 50% of the energy flow in some coral reef habitats (Herler 
et al., 2011). Despite their evolutionary and ecological impor- 
tance, the phylogenetic relationships among species within 
Gobiidae and their location within Gobioidei are still poorly 
understood (Murdy, 1989; Parenti and Thomas, 1998; Thack- 


er and Schaefer, 2000; Larson, 2001). To date, the classifi- 
cation of gobies still remains largely reliant on external mor- 
phology (Pezold, 1993; Akihito et al., 2000; Nelson, 2006), 
and diagnostic characters separating species are subtle and 
problematic. 

Molecular biology has contributed to addressing taxon 
identification and phylogenetic relationship questions. Mito- 
chondrial DNA (mtDNA) markers have historically formed 
the core of most molecular systematic analyses and are still 
the most widely used for reconstructing phylogeny (Brown 
et al., 1979; Moore, 1995; Johns and Avise, 1998); this is 
probably due to their single copy nature and relative ease of 
sequencing (Moore, 1995). Genetic divergence is also enhanc- 
ed by the higher rate of sequence evolution in vertebrate 
mtDNA compared to that of nuclear coding regions (Johns 
and Avise, 1998). However, the choice of a suitable gene is 
crucial for identification and phylogenetic reconstruction 
among closely related species (Brown et al., 1979; Moore, 
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1995; Johns and Avise, 1998), because different parts of the 
mtDNA genome evolve at different rates (Avise and Ellis, 
1986; Roques et al., 2006). 

Mitochondrial cytochrome oxidase subunit I (COI) could 
serve as a rapid and reliable barcoding marker for identify- 
ing species and for discovering new species across the entire 
animal kingdom (Hebert et al., 2003). Although skepticism 
has frequently been expressed (Ebach and Holdrege, 2005; 
Will et al., 2005), DNA barcoding based on COI has been 
successful to identify species across a wide array of taxa over 
the last decade (e.g., Hebert et al., 2004; Clare et al., 2007; 
Hubert et al., 2008; Feng et al., 2011). A clear gap should 
exist between intra- and interspecific COI sequence diver- 
gence with about a 20-fold difference for DNA barcoding to 
be perfectly effective in delimitating species (Hebert et al., 
2003). A standard sequence threshold can be projected to 
outline species boundaries by employing this barcoding gap. 
However, utilizing such a threshold value may be challeng- 
ing, particularly when attempts do not include numerous spe- 
cimens, such as for critically endangered taxa. 

Only a few studies have addressed gobioid interrelation- 
ships based on molecular data (e.g., Akihito et al., 2000; 
Wang et al., 2001; Thacker, 2003, 2009; Thacker and Hard- 
man, 2005). These studies used different taxon and nucleo- 
tide sampling methods. Yet, testing the effectiveness of COI 
DNA barcodes on species that have emerged as a result of 
radiation, such as gobies, remains a major challenge in evo- 
lutionary biology. Here, we sequenced the COI of 48 spe- 
cies collected from South Korea (Table 1) to evaluate the 
performance of distance-based DNA barcoding for phylo- 
genetic analyses. We specifically aimed to provide novel 
data on a comparison of pairwise divergence levels among 
species in the same genus vs. species in different genera. 
GenBank sequences were also included in the analyses to 
use a dataset with large taxonomic coverage (n “144 species) 
(Table 2). 


MATERIALS AND METHODS 


Sample collection 

Fish were collected using seine and dip nets from January to 
November 2011 from 21 sites across freshwater systems, 
coastal areas, and the ocean near South Korea (Table 1). Spe- 
cimens were identified based on morphological characters. 
Entire bodies of all individuals were preserved in 95% etha- 
nol, and 44 nominal species were sequenced for COI gene 
fragments. 


DNA isolation, amplification, and sequencing 
We used the Wizard Genomic DNA purification kit (Promega, 
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Madison, WI, USA) to extract genomic DNA from the right 
pectoral fin of each fish specimen. The COI was amplified 
using gobiid-specific primers: GOBYF7558 (forward) 5’- 
TIT GCW ATT ATG GCW GGA TTT G-3' and GOBYB 
8197 (reverse) 5’-ATT ATT AGG GCG TGG TCG TGG-3’ 
(Thacker, 2003) and COI fish universal primers, FF2d (for- 
ward) 5'- TTC TCC ACC AAC CAC AAR GAY ATY GG- 
3’ and FF1d (reverse) 5'-CAC CTC AGG GTG TCC GAA 
RAA YCA RAA-3’ (Ivanova et al., 2007). Each polymerase 
chain reaction (PCR) amplification was carried out in a 50 
uL reaction volume composed of ~75 ng DNA extract, 0.25 
mM of each deoxynucleotide, 0.25 mM of each forward and 
reverse primer, 3mM MgCh, 1 x PCR buffer, and 0.25 units 
of Taq DNA polymerase (Solgent, Daejeon, Korea). GenePro 
(BIOER) was used to amplify the COI with the following 
program: 94°C for 10 min, 35 cycles of 30s at 94°C, 30 s at 
54°C (for GOBYF7558 - GOBYB8197) and 52°C (for FF2d- 
FF1d), 30s at 72°C and final elongation at 72°C for 10 min. 
PCR products were loaded on 1% agarose gels containing 
0.003% ethidium bromide and visualized using the GelDoc- 
It TM Imaging System (UVP). Amplifications were consid- 
ered successful when a expected sized band was observed 
on the agarose gel. PCR products were cleaned using a PCR 
purification kit (Solgent). The COI was sequenced directly 
using the BigDye-Terminator V3.1 kit (Applied Biosystems, 
Foster City, CA, USA) and an ABI3730XL sequencer at 
Genotech (Daejeon, Korea). 


Sequence data analyses 

Complementary DNA sequences were assembled using the 
Bioedit 5.0.9 sequence-editing software (Hall, 1999). Sequ- 
ences were aligned using Clustal X 2.0 default settings (Lar- 
kin et al., 2007). Alignments were translated to amino acids 
under the vertebrate mitochondrial option using MEGA 5 
(Tamura et al., 2011) to detect frameshift mutations and pre- 
mature stop codons, which may indicate the presence of 
pseudogenes. The Genbank accession numbers of newly 
determined sequences were JX679021 -JX679066 and are 
listed in Table 1. Genetic distances were calculated to quan- 
tify sequence divergences among species using both p dis- 
tance and Kimura two-parameter (K2P) models (1,000 boot- 
strapping) (Kimura, 1980), as implemented in MEGA. Rates 
of synonymous and nonsynonymous substitutions were also 
calculated with MEGA using both standard and modified (at 
1.4 standard errors with 1,000 bootstrapping samples) Nei- 
Gojobori models (Nei and Gojobori, 1986; Nei and Kumar, 
2000). Genetic distances were calculated at intrageneric, 
intrasubfamilial, intrafamilial, and interfamilial levels. Alto- 
gether, 10,296 pairwise distances were compared in this study. 
The degree of sequence conservation per site, Ryeg, was 
defined as Ryeg=2—(— 22plogəp) (Ward and Holmes, 2007), 
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Table 2. List of reference species used in this study 


(Sub)family Species Reference Gadb 
accession no. 
Odontobutidae Perccottus glenii Thacker and Hardman (2005) AY722171 
Rhyacichthyidae Protogobius attiti Keith et al. (2011) HQ639032 
Rhyacichthys guilberti Keith et al. (2011) HQ639030 
Benthophilinae Neogobius melanostomus Unpublished HQ960511 
Proterorhinus marmoratus Hubert et al. (2012) EU524305 
Proterorhinus semilunaris Unpublished HQ961006 
Gobiinae Amblyeleotris aurora Unpublished HQ561508 
Amblyeleotris guttata Steinke et al. (2009) FJ582710 
Amblyeleotris steinitzi Steinke et al. (2009) FJ582712 
Amblyeleotris wheeleri Unpublished HQ561506 
Amblygobius decussatus Steinke et al. (2009) FJ582722 
Amblygobius phalaena Hubert et al. (2012) JQ431409 
Asterropteryx ensifera Hubert et al. (2012) JQ431470 
Asterropteryx semipunctata Hubert et al. (2012) JQ431471 
Bathygobius coalitus Hubert et al. (2012) JQ431479 
Bathygobius cocosensis Hubert et al. (2012) JQ431480 
Bathygobius cotticeps Hubert et al. (2012) JQ431482 
Bathygobius curacao Weigt et al. (2012) JQ839961 
Bathygobius lacertus Weigt et al. (2012) JQ839962 
Bathygobius laddi Unpublished JF492946 
Bathygobius mystacium Weigt et al. (2012) JQ839963 
Bathygobius soporator Tornabene et al. (2010) HM748425 
Caffrogobius caffer Unpublished HQ945885 
Coryphopterus tortugae Weigt et al. (2012) JQ842827 
Croilia mossambica Unpublished HQ561481 
Cryptocentrus cryptocentrus Unpublished HQ561467 
Cryptocentrus leptocephalus Steinke et al. (2009) FJ583296 
Cryptocentrus pavoninoides Steinke et al. (2009) FJ583298 
Elacatinus evelynae Steinke et al. (2009) FJ583388 
Elacatinus oceanops Steinke et al. (2009) FJ583389 
Eviota afelei Thacker (2003) AF391391 
Eviota disrupta Leray et al. (2012) JN107907 
Eviota distigma Hubert et al. (2012) JQ349971 
Eviota indica Hubert et al. (2012) JQ349972 
Eviota prasina Hubert et al. (2012) JQ349973 
Favonigobius exquisitus Thacker et al. (2011) HQ909465 
Fusigobius signipinnis Steinke et al. (2009) FJ583414 
Glossogobius aureus Aquino et al. (2011) HQ682689 
Glossogobius callidus Steinke et al. (2009) JF493535 
Gobiodon ceramensis Steinke et al. (2009) FJ583428 
Gobiodon histrio Steinke et al. (2009) FJ583450 
Gobiodon okinawae Steinke et al. (2009) FJ583454 
Gobiodon quinquestrigatus Hubert et al. (2012) JQ431768 
Gobiopterus lacustris Aquino et al. (2011) HQ682693 
Gobius bucchichi Unpublished JF935258 
Gobius cruentatus Unpublished JF935263 
Istigobius decoratus Steinke et al. (2009) JF493692 
Istigobius rigilius Thacker et al. (2011) HQ536672 
Knipowitschia caucasica Triantafyllidis et al. (2011) HQ600736 
Nes longus Weigt et al. (2012) JQ841296 
Paragobiodon lacunicolus Steinke et al. (2009) FJ583828 
Pomatoschistus tortonesei Unpublished FJ751922 
Priolepis cincta Hubert et al. (2012) JQ350251 
Priolepis compita Hubert et al. (2012) JQ432034 
Priolepis eugenius Thacker (2003) AF391329 
Priolepis hipoliti Thacker et al. (2011) HQ909484 
Priolepis inhaca Hubert et al. (2012) JQ432036 
Priolepis semidoliata Hubert et al. (2012) JQ432038 
Priolepis squamogena Hubert et al. (2012) JQ432040 
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(Sub)family Species Reference GonBanle 
accession no. 
Gobiinae Priolepis triops Hubert et al. (2012) JQ432042 
Pycnomma roosevelti Unpublished GU224480 
Rhinogobiops nicholsii Unpublished JN582118 
Trimma caesiura Thacker (2009) EU381039 
Trimma macrophthalma Hubert et al. (2012) JQ350401 
Trimma mendelssohni Hubert et al. (2012) JQ350407 
Trimma milta Hubert et al. (2012) JQ432200 
Valenciennea helsdingenii Steinke et al. (2009) FJ584201 
Valenciennea longipinnis Steinke et al. (2009) FJ584202 
Valenciennea muralis Steinke et al. (2009) FJ584203 
Valenciennea puellaris Steinke et al. (2009) FJ584212 
Valenciennea sexguttata Steinke et al. (2009) FJ584224 
Valenciennea strigata Unpublished HQ945877 
Valenciennea wardii Steinke et al. (2009) FJ584240 
Gobionellinae Awaous aeneofuscus Unpublished HQ945950 
Awaous melanocephalus Aquilino et al. (2011) HQ654674 
Ctenogobiops tangaroai Steinke et al. (2009) FJ583310 
Ctenogobius saepepallens Weigt et al. (2012) JQ840025 
Gnatholepis cauerensis Unpublished HQ561511 
Gnatholepis thompsoni Weigt et al. (2012) JQ840079 
Gobioides broussonnetii April et al. (2011) JN026727 
Gobionellus oceanicus Weigt et al. (2012) JQ841902 
Lethops connectens Unpublished GU440372 
Oligolepis acutipennis Aquilino et al. (2011) HQ654730 
Oligolepis keiensis Unpublished HQ945926 
Stenogobius polyzona Unpublished HQ945939 
Typhlogobius californiensis Unpublished GU440562 
Oxudercinae Scartelaos histophorus Unpublished FJ238032 
Taenioides sp. Steinke et al. (2009) FJ584167 
Sicydiinae Akihito vanuatu Keith et al. (2011) HQ639065 
Cotylopus rubiripinnis Keith et al. (2011) HQ639038 
Lentipes armatus Keith et al. (2011) HQ639070 
Sicydium punctatum Keith et al. (2011) HQ639050 
Sicyopterus lagocephalus Hubert et al. (2012) JQ432152 
Sicyopterus pugnans Hubert et al. (2012) JQ432154 
Sicyopus chloe Keith et al. (2011) HQ639058 
Stiphodon elegans Hubert et al. (2012) JQ432172 


where p is the observed frequency of each base at a particu- 
lar position and the maximal degree of conservation was 2, 
which was achieved when all nucleotides at a particular site 
in the 144 species were the same. 

A Bayesian inference (BI) tree was established using 
MrBayes 3.1.2 (Ronquist and Huelsenbeck, 2003) with two 
outgroup species from the Family Rhyacichthydae, Proto- 
gobius attiti and Ryacichthys guilberti. The best-fit model 
of DNA sequence evolution was chosen using ModelTest 
3.8 (Posada and Crandall, 1998) and Akaike information 
criteria; the chosen model was GTR+I+G. The analysis 
was run for 10 million generations with sampling of one tree 
every 500 generations. Two independent Markov Chain 
Monte Carlo runs were conducted simultaneously. The first 
1,000 trees of each run were discarded as burn-in. 
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RESULTS 


The COI genes of each species were confidently aligned, 
and the equivocal bases at each end were trimmed to yield a 
final sequence of 542 bp. No indels were detected. Transla- 
tion of the sequences did not reveal frame-shift mutations or 
premature stop codons, confirming that our amplified frag- 
ments were functional. Among the 542 nucleotide positions, 
245 were polymorphic, and 230 were parsimony informative. 
The proportion of T, C, A, and G bases for all 144 sequences 
was 30.596, 28.196, 23.196, and 18.296, respectively. The GC 
content was relatively higher at the first codon base (56.396) 
than that at the second (43.1%) or third (39.6%). The degree 
of conservation (y.q) was calculated for each base of the 542 
nucleotides; the most common and maximum value was 2, 
which was achieved when all nucleotides at a particular site 
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Table 3. Estimated evolutionary parameters ( x 100) for the nucleotide substitutions in the cytochrome oxidase c subunit 1 (COI) 


barcoding region from 144 gobioid fish species 


All fish Gobiinae Gobionellinae 
Genetic distance 
Mean SD Mean SD Mean SD 
Codon nucleotide position 
p distance First 5.51 1.63 5.45 2.05 5.21 2.06 
Second 0.18 0.15 0.28 0.47 0.02 0.11 
Third 54.31 3.61 55.49 5.19 50.19 7.05 
Overall 19.94 1.66 20.34 2.10 18.42 2.69 
K2P distance First 6.23 2.09 6.14 2.53 5.85 2.52 
Second 0.18 0.15 0.28 0.48 0.02 0.11 
Third 403.10 439.15 499.53 552.80 240.11 145.15 
Overall 28.18 4.06 29.06 4.05 25.34 4.42 
Synonymous and nonsynonymous distances 
Nei-Gojobori (N-G) Synonymous 71.91 6.98 73.41 6.77 66.95 9.32 
Nonsynonymous 1.64 0.94 1.72 0.99 1.29 0.83 
dw/ds 0.02 0.00 0.02 0.00 0.02 0.00 
Modified N-G Synonymous 62.78 6.06 64.06 5.88 58.44 8.11 
Nonsynonymous 1.73 1.00 1.81 1.05 1.36 0.88 
dw/ds 0.03 0.00 0.03 0.00 0.02 0.00 


K2P, Kimura-2-parameter. 


of the 144 species were the same. Every third codon base 
was highly variable with a 0.64 mean R,4. The first (1.72) 
and the second codon bases (1.98) were nearly monomorphic. 
Nucleotide genetic distance parameters, p and K2P distance, 
also showed an almost zero rate of substitution for second 
nucleotide positions, with the first position being an order 
of magnitude higher and the third position being unparallel- 
ed among them (Table 3). The rate of synonymous substitu- 
tions was much higher than the rate of nonsynonymous sub- 
stitutions (Table 3). 

Three species could not be separated using the COI sequ- 
ence analysis, including Chaenogobius gulosus (Gobionelli- 
nae), Chaeturichthys stigmatias (Gobionellinae), and Lophio- 
gobius ocellicauda (Gobiinae). The C. gulosus sequence was 
highly similar to those of other congeneric species, as expect- 
ed, whereas the C. stigmatias and L. ocellicauda results were 
very surprising. Multiple specimens of these species should 
be extensively analyzed in a future study to check the genet- 
ic divergence among these species; thus, those sequences 
were not deposited in GenBank. The average diversity among 
142 haplotypes (Ha) was 0.999 + 0.001 (mean + standard de- 
viation), and the average nucleotide diversity (7) was 0.199 
2: 0.002. Twenty-six genera were represented by two or more 
species. Levels of intragenus divergence were generally high 
(Table 4) but varied greatly among genera. For example, the 
average within-genus divergence of Awaous, Eviota, and 
Trimma was 29.85%, 29.13%, and 28.71%, respectively, 
which was larger than that of overall interspecific divergence. 
These values were considerably larger than those of Amblygo- 
bius, Asterropteryx and Rhinogobius (9.76%, 9.87%, and 
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Table 4. Mean genetic divergences for the cytochrome oxidase 
c subunit 1 (COI) nucleotide sequences (p and Kimura-2-para- 
meter [K2P] distances) among 144 gobioid species 


p distance ( x 100) K2P distance ( x 100) 


Taxonomic level 


Mean Min Max Mean Min Max 
Within-genus 15.63 0.18 23.80 20.87 0.19 36.13 
Within-subfamily 19.69 9.20 26.20 27.87 0.19 41.44 
Within-family 20.05 6.83 26.94 28.35 0.94 44.53 
Among-family 20.18 10.89 25.83 28.56 13.27 41.06 


Four taxonomic levels are represented such as within-genus, within-sub- 
family, within-family, and among-family. 


10.74%, respectively). Mean interspecific divergences at 
higher taxonomic levels were slightly larger than that at the 
intragenus level, resulting in large overlaps among levels 
(Table 4). As roughly 90% of all gobiid fishes are either in 
Gobiinae or Gobionellinae, these two subfamilies were com- 
pared for nucleotide substitution rate (Table 3). Species from 
Gobiinae were consistently higher than Gobionelline fishes 
in every parameter estimated (Table 3). 

Our phylogenetic data provide little evidence to support 
the previous claims at the generic and higher taxonomic 
levels, based on phenetic analyses. Several species did not 
cluster into their respective groups, and the BI tree failed to 
correctly identify some genera or subfamilies (Fig. 1). For 
example, Tridentiger barbatus (Gobiinellinae) clustered into 
the Asterropteryx clade (Gobiinae) with high nodal support 
rather than with other species of Tridentiger. Pterogobius 
zacalles clustered with Acanthogobius, whereas P. zonoleu- 
cus formed a monophyletic group with the Rhinogobius brun- 
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Fig. 1. Molecular phylogeny of Gobioidei based on Bayesian inference from 542 bp of the mitochondrial cytochrome oxidase c sub- 
unit 1 (COI) gene with two outgroup species from the Family Rhyacichthyidae: Protogobius attiti and Rhyacichthys guilberti. Under 
the chosen model, GTR+I+G, analysis was run for 10 million generations with sampling of one tree every 500 generations. Numbers 
above branches indicate posterior probabilities (7 0.8). Abbreviations in brackets on the right side indicate higher taxonomic names 
(subfamily and family) in current usage, including Amblyopinae (A), Benthophilinae (B), Gobiinae (I), Gobionellinae (N), Oxudercinae 
(O), Sicydiinae (S), Microdesmidae (MD), Odontobutidae (OD), and Rhyacichthyidae (RD). Asterisks immediately after the higher 
taxonomic names and arrows on nodes indicate taxa failing to resolve monophyletic assemblages, and species clustered into unre- 
lated groups, respectively. 
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neus complex, albeit without good nodal support. The BI tree 
also failed to resolve monophyletic assemblages of some 
taxa, such as Awaous, Eviota, Trimma, and Ondontobutidae. 
Those exceptions aside, the tree largely assigned species to 
identical major groups. 


DISCUSSION 


The COI sequences did resolve genealogical relationships 
well at the level of genera and family in Gobioidei. As pre- 
viously noted in Che et al. (2012), this was possibly due to 
at least two factors. Most likely, the 542 nucleotides do not 
provide sufficient phylogenetically informative characters 
to recover the true phylogeny when examining hundreds of 
taxa with enormous diversity. In addition, the fast mutation 
rate and saturation in the third codon position can be a dis- 
advantage at deeper phylogenetic levels, and the subsequent 
long terminal branches may impede resolution of ancient 
speciation due to the chance accumulation of shared charac- 
ter states (Huelsenbeck, 1997). Despite the poor monophylet- 
ic resolution in several taxa, some clear phylogenetic signa- 
tures were observed in the COI sequence data. For example, 
several major congeneric species, Rhinogobius and Bathygo- 
bius, tended to cluster together with no exception and, in 
most cases, so did consubfamilial species. It is believed that 
the utilization of more gene regions including nuclear DNA 
will assist in offering a more reliable phylogeny within Gob- 
iidae and their placement within Gobioidei. Several nuclear 
genes such as recombination activating genes 1 and 2 (Rag 1 
and 2) and ryandine receptor 3 (Ryr3) may have slower rate 
of sequence evolution in gobies compared to that of mtDNA 
genes (Yamada et al., 2009). 

Our results support the validity of COI barcoding for spe- 
cies identification in gobiid species, although no attempt 
was made to include numerous specimens for any one spe- 
cies. One fundamental barcoding criterion is that congeneric 
divergence should be significantly higher than that of con- 
specific divergence (Hubert et al., 2008). The average intra- 
genus distance (K2P) for 28 genera with multiple species in 
the present study was 21.09%, which was considerably high- 
er than the values obtained among fish species in previous 
studies (9.93% from Ward et al., 2005; 9.54% from Ward 
and Holmes, 2007). In addition, the interspecific divergence 
at higher taxonomic levels was not significantly larger than 
that at the intragenus level, suggesting that congeneric gobies 
possess substantial interspecific sequence divergence in their 
COI genes. Significantly more nucleotide changes were 
observed at the third codon position than those at the first 
and the second, revealing that the extensive variation shown 
among the COI sequences typically reflects synonymous 
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changes and little variation at the protein level. Consequent- 
ly, the proportion of nonsynonymous to synonymous changes 
was far less than one (Table 2). As previously noted in Ward 
and Holmes (2007), this result must be due to exceptionally 
strong purifying selection of the COI gene and confirms that 
the ability of COI to identify species in Gobiidae is depen- 
dent on the degenerate nature of the genetic code. 

We did not provide sufficient evidence for the utility of the 
COI towards cryptic species identification in several species 
complexes. Gobiidae taxonomy has been studied extensive- 
ly for the last several decades, but confusion still exists. One 
typical case is the taxonomic status of the Gymnogobius spe- 
cies complex (G. urotaenia, G. opperiens, and G. petschilien- 
sis) (e.g., Harada et al., 2002). In our results, overall K2P 
divergence within Gymnogobius was 12.7496, whereas the 
average value among the Gymnogobius species complex 
was only 1.81%, probably reflecting a short history of repro- 
ductive isolation. The Rhinogobius brunneus complex (in 
our analysis, Rhinogobius brunneus, R. brunneus CB, and 
R. brunneus CO) is also a representative example of a cryp- 
tic species complex in Gobiidae (Kawanabe and Mizuno, 
1989; Kim, 1995). Overall, K2P divergence within Rhino- 
gobius was > 10%, whereas the average value among the R. 
brunneus complex was just 226. Although more work needs 
to be done with multiple specimens, the COI sequence may 
not be a reliable tool to delineate cryptic and complex species 
boundaries in the family Gobiidae. 
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